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Abstract 

The muon and anti-muon neutrino energy spectrum is determined from 2000- 
2003 AMANDA telescope data using regularised unfolding. This is the first 
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measurement of atmospheric neutrinos in the energy range 2 - 200 TeV. 
The result is compared to different atmospheric neutrino models and it is 
compatible with the atmospheric neutrinos from pion and kaon decays. No 
significant contribution from charm hadron decays or extraterrestrial neutri- 
nos is detected. The capabilities to improve the measurement of the neutrino 
spectrum with the successor experiment IceCube are discussed. 

Key words: atmospheric neutrinos; unfolding; neural net; AMANDA; 
Cherenkov radiation 

PACS: 95.55.Vj, 95.85.Ry, 29.40.Ka, 14.60.Lm 



1. Introduction 



At energies above 0.1 TeV, about one cosmic ray particle per square meter 
per second reaches Earth. At the highest observed energies, particles reach 
more than 10^° eV, which is far above what can be achieved in man-made 
accelerators. The origin of these charged cosmic rays is still being discussed, 
as their direction is scrambled by extragalactic and galactic magnetic fields. 
One option for identifying the origin of cosmic rays is the observation of 
secondary particles produced in cosmic ray interactions in the astrophysical 
plasmas themselves: if a proton p interacts with ambient matter or photon 
fields 7, pionic secondaries are produced vi a the processes p + p — ?■ vr + X and 



p+7 — 7- A"^ — 7- n + 7r''"/p+7r°, respectively (lAmsler et al.l . l2008l ). The charged 



pions subsequently decay into neutrinos, vr^ — )■ jji^ 



where we do not distinguish between neutrinos and anti-neutrinos. The re- 
sulting neutrino fiux usually follows the spectral behaviour of the protons, 
whic h is predicted to be close to dN/ dE oc according to Fermi accelera- 
tion fiFermil f ll949l . Il954[ )l The conventional atmospheric neutrino spectrum 



approximately dN/dE oc E ^■'^ { 


Honda et al. . 


1995 


Volkova and Zatsepin, 


1980; 


Gaisser et al.. 


2001; 


Barr et al.. 


2004; 


Honda et al.. 


2004 


). An addi- 



tional component of the atmospheric neutrino fiux comes from the decays 
of hadrons containing charm and bottom quarks. This fiux, known as the 
prompt component is expected to have a spectrum close to dN/ dE oc 



fiBugaev et al.l.ll989l : ICostal . l200ll : iFiorentini et all . l200ll : Martin et al.l . 12003 



Honda et al.l . |2004| . e.g.). The prompt atmospheric neutrino fiux is lower 



than the conventional fiux but could start to dominate the total spectrum 
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at energies a bove about 100 TeV. So f ar, only the conventional neutrino flux 



is observed (jAchterberg et al.l . l2007al ). Measurements at high energies, i.e. 



above 10 — 100 TeV, provide an opportunity to reveal an extraterrestrial or 
a charm component. At these energies, the Antarctic Muon And Neutrino 
Detector Array (AMANDA) and its successor IceCube are able to make mea- 
surements to look for deviations from the conventional atmospheric neutrino 
flux. 

The AMANDA-II detector was designed for the detection of neutrinos above 
100 GeV. It is composed of 677 Optical Modules (OMs), each containing a 
8-inch, 14-dynode photomultiplier tube (PMT) and a voltage divider for the 
high voltage. The PMTs are optically coupled to the pressure glass sphere 
with a silicon gel and can be operated at a high gain of about 1 ■ 10^. The 
optical modules are attached to 19 vertical strings, instrumenting a cylindri- 
cal volume of 0.016 km^ ( with a radius of 100 m and a height of 500 m), see 
e.g. iDe Young et al.l (12008 ). Secondary muons in the ice are produced via the 
process U/j^ + N ^ fi + Xu- The muons produce Cherenkov radiation if they 
travel faster than the speed of light in ice (i.e. if the muons travel faster than 
V > 0.8 ■ Co, with Co as the speed of light in vacuum). Additional Cherenkov 
radiation comes from the particles produced in muon interactions, such as 
bremsstrahlung, direct pair production and photonuclear interactions, all 
dominating at muon energies above 1 TeV. At higher energies, the sum of 
the energy loss due to stochastic processes (i.e. bremsstrahlung, pair pro- 
duction and nuclear interaction) is dominant and increases linearly with the 
energy. The amount of light detected with the optical modules rises with 
the muon energy and therefore also with the energy of the parent neutrino. 
Thus the detected light amount can be used to determine the primary neu- 
trino energy spectrum. Neutrino-induced muons can be distinguished from 
atmospheric muons by selecting events that traverse the Earth and arrive at 
the detector from below the horizon. Atmospheric muons cannot reach the 
detector from those directions since they are absorbed on their way through 
the Earth. In this energy range {E < 200 TeV), neutrino absorption in the 
Earth is not signiflcant. Neutrinos can traverse the matter without loss and 
some neutrinos interact close to the detector, so that the products of these 
neutrino interactions can be observed. 



^Here and throughout the paper, we use the same notation for particles and antiparti- 
cles. 
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AMANDA data from the years 2000 to 2003 are analyzed to determine the 
energy spectrum of neutrinos, presenting for the first time the atmospheric 
neutrino spectrum in the energy range 2 — 200 TeV. In section [21 predictions 
for atmospheric neutrinos are reviewed. In section [3l a conceptual overview 
of the issues involved in deconvolving a spectrum from observed data are dis- 
cussed. In section HI more details of the data reduction, simulation and anal- 
ysis method for the deconvolution of the neutrino spectrum are explained. 
Section [5] then describes a neural network used for the construction of an op- 
timal energy-correlated variable, while section |6] shows how the atmospheric 
spectrum is determined by regularised unfolding, and discusses the sources 
of statistical and systematic uncertainties that enter the calculation. Section 
[7] summarises the results while section |H] discusses them in the context of 
other experimental results and flux predictions. Finally, section [9] gives the 
conclusions from this analysis along with an outlook on the possibilities for 
IceCube, the successor of the AMANDA experiment. 



2. Atmospheric neutrinos 



When cosmic rays traverse the Earth's atmosphere, hadronic showers are 
produced by their interactions with the atmosphere. Depending on the en- 
ergy of the primary cosmic ray, different secondaries can be produced. Up 
to energies of ~ 100 TeV, the flux is dominated by pion and kaon decays. 
This flux is usually referred to as the conventional atmospheric neutrino flux. 
The expected power-law spectrum of conventional atmospheric neutrinos is 
typically one power steeper than the primary cosmic ray spectrum: charged 
pions and kaons have rest-frame lifetimes of the order of lO"*^ s. High energy 
pions and kaons travel far enough in the atmosphere that they may interact 
with atmospheric nuclei before they decay. As the lifetime increases with the 
particles' energy, high-energy pions and kaons have a higher probability to 
interact before decaying, which steepens the spectrum by one power. An ana- 
lytic descri ption of the neutrino spectrum between 100 GeV and 5.4-10^ GeV 
is given by IVolkova and ZatsepinI (Il980l ) : 



dN 



dE^dn 



0.213 



l.UE,/EK±{e)\ 



(1) 



with = 0.0285 GeV"^cm-2s-^sr-^ and 7 = 2.69. Here, E^ a nd Ek± are 

energy distribution parameters that depend on the zenith angle 6 fIVolkova and Zatsepin 
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1980|). Equation ([T]) covers the energy range relevant for this analysis, and 
will be used to fit the observed spectrum. 

The energy spectrum varies with the zenith angle, as a pion traveling through 
the atmosphere horizontally (cos(^) = 1) experiences a smaller density gra- 
dient than a pion traversing the atmosphere vertically (cos(6') = 0). Thus, 
nearly vertical pions have a higher probability of interacting with the at- 
mosphere, which reduces the flux compared to the horizontal component. 
Other predict i ons o f conventiona l atra osphe ric neutrino f l uxes are giv en by 
Gaisser et all fl200lh : iBarr et al.l f l2004h and iHonda et all jgOOj . l2007h . with 



uncertainties in the modeling of around 15% ( Barr et al. . 20061 ) 



The prompt atmospheric neutrino flux is a second component in the atmo- 
spheric neutrino spectrum and is due to the decay of charm and bottom 
hadrons, which contain charm quarks. Since D~ and hadrons have life- 
times shorter than 10~^^ s, they decay before any further interaction with 
the ambient matter can take place. Thus, these hadrons produce an isotropic 
neutrino spectrum whose shape is close to the primary cosmic ray spectrum, 
i.e. E- 



-2.7 



The atmospheric neutrino energy spectrum has so far been measured in the 



range 1—10 TeV ( iDaum et al.l . ll995uGonzalez-Garcia et al.l . l2006l : lAbbasi et al 



20091 ). The flux is found to follow the prediction of the conventional atmo- 
spheric flux within uncertainties. Here, the atmospheric neutrino spectrum 
extending up to 200 TeV is presented for the first time. As theoretical uncer- 
tainties increase towards higher energies, the experimental investigation of 
the atmospheric spectrum can lead to conclusions about particle interactions, 
possible charm contribution or a possible extraterrestrial component. 



3. Determination of the neutrino energy spectrum 

Determination of the neutrino energy spectrum in a detector like AMANDA 
is complicated by various factors. Neutrino energies are not measured directly 
but are inferred from measuring the energies of the interaction products of 
the neutrinos. The three flavors of neutrinos produce different secondary 
patterns of particles in a detector, with correspondingly different correlations 
to the primary energy. An electron-neutrino will deposit most of its energy 
into an electromagnetic or hadronic cascade, which, limited by the resolution 
of the energy reconstruction algorithm, will correlate directly to the primary 
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energy. A tau-neutrino produces two cascades of particles, one from the 
initial interaction and a second from the decay of the tau particle. If the 
energy of the event is sufficiently small, then both of these may be contained 
in the detector and thus be representative of the neutrino energy. For the 
main signature in AMANDA, upward moving muon tracks, the neutrino 
energy cannot be directly measured, because of the range of the muon. For 
a neutrino that interacts inside the detector, the muon may carry away a 
significant fraction of the neutrino energy which is then not measured. For 
high-energy neutrinos that interact very far away, the energy deposited in 
the detector from the final muon will only correlate weakly with the neutrino 
energy. This further degrades the correlation to the initial energy. Moreover, 
whatever energy algorithm is used, its limited resolution will further degrade 
the ability of the experiment to measure a particle spectrum. 

In practise, one uses knowledge of the neutrino physics and detector response 
to infer the parameters of the primary neutrino spectrum. Mathematically, 
this inference is represented by an integral equation. 



where p{y) is the probability of observing an event with a reconstructed pa- 
rameter vector y, A{y, E) represents the detector response to initial particles 
of energy E, is the neutrino fiux and b{y) is the non-neutrino back- 

ground. Here, p{y) represents the distribution of some parameter related 
to the energy of the event, e.g. the total deposited energy, or a direct esti- 
mate of the particle energy. Once p{y) is obtained and A{y, E) is determined 
via detector calibration and simulation, then there are several ways to infer 
$(£'). The simplest is to take theoretical estimates of and see which 

one results in an expected p{y) that fits the data best. If the theory can be 
parameterised, then these parameters can be adjusted until a best fit is found 
and errors on the parameters can be determined. Examples of this would be 
direct fitting of the slope 7 and normalisation C of a spectrum of the form 
^{E) ~ C ■ E-^. 

Another method employed in AMANDA is to assume the general form of a 
theoretical fiux from a calculation and fit for a free normalisation and de- 
viation from the spectral shape. These methods are commonly known as 
forward folding - the parameters of the fiux are adjusted and then forward 
propagated through to an observable which is compared with data. In prin- 
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ciple, the flux can be made theory independent by choosing a parametric 
form with a large number of parameters. This leads to a problem where 
small statistical fluctuations in the observed data may result in unphysical 
solutions for the spectrum, manifesting as bumps and dips in the spectrum. 
Methods known as regularisation are employed, whereby constraints on the 
smoothness of the solution are imposed, to control these effects. For the low 
parameter methods, regularisation is built into the solution, for instance by 
the assumption of a power law spectral form. 

As the number of parameters grows, the computational time to solve for 
the best flt set of parameters also increases. As an alternative, non-iterative 
methods, collectively known as unfolding methods, exist for the direct so- 
lution of the high parameter problems. The difficulties of regularisation 
are still inherent in these direct methods, just as they are in the intera- 
tive forward folding methods. These considerations aside, the conceptual 
basis of unfolding is to use the direct inverse of the matrix A{y, E) to solve 
for$(E) = A(2/,E)-V(y). 

The accuracy of a deconvolution of the spectrum using any of the meth- 
ods improves on finding an observable y (which is possibly vector-valued) 
that is well correlated to the neutrino energy, manifested by minimising the 
infiuence of the off-diagonal elements of the matrix A{y, E). There is an irre- 
ducible component of the off-diagonal elements set by the physics of neutrino 
interaction and muon propagation. For instance, the fact that the muon car- 
ries away only a fraction of the neutrino energy and then loses energy as it 
propagates to the detector reduces the correlation of the muon energy to the 
original neutrino energy. The size of the off-diagonal elements depends on 
how well correlated the chosen variables are with the neutrino energy. In this 
work, the energies of the muons are reconstructed by using a neural network, 
making use of six energy sensitive variables. The output of the neural net- 
work, along with two other energy sensitive parameters, form the observable 
y for the unfolding. 

We do not directly unfold the atmospheric neutrino spectrum but unfold 
an intermediate spectrum which is a convolution of the actual flux with 
the probability of the neutrinos having passed through the Earth and with 
the efficiency for detection. Thus, the intial unfolding returns a neutrino 
spectrum with the same number of events as observed in the detector. This is 
then corrected as a function of energy and angle back to the true atmospheric 
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spectrum. 



4. The data set and the analysis method 



For the determination of the neutrino energy spectrum from 2 — 200 TeV, 
807 days effective hvetime of data taken by the AMANDA detector between 
the years 2000 and 2003 are used. The selection of neutrino candidates for 
zenith angles 9 > 90° is presented in lAchterberg et al.l ( l2007al ). In a final step, 
tracks with 6 < 100° are removed to minimise the atmospheric muon con- 
tamination of t he neutrino sample. The final sample contains 2972 neutrinos 
( iMiinichl . 120071 ) and includes a background of less than 1% misreconstructed 
atmospheric muons ( lAchterberg et al.l . |2007aJ) . 



Neutrino 
generator 
nusim 



Air stiower 

generator 

CORSIKA 



Muon propagation 
MMC 



Detector MC 
AMASIM 



Energy 
variables 



CH 
CHI 



HITS 



RMS 



Neural Net 



Unfolding 



Additional 
Variables 



log(N(,j^) 

log( Q ) 



Event Spectrum 



Particle 
Spectrum 



Figure 1: Scheme of the simulation and analysis chain 



In order to avoid possible biases we perform a blind analysis in the following 
sense: The properties of the selected events are initially checked on 10% of 
the full data set, by comparing this reduced data set to simulations of the 
atmospheric muon background and the conventional atmospheric neutrino 
flux. The simulation is done by using the air shower simulation CORSIKA 
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for the background (lHeckl.ll998r) an d the neutrino generator nusim for the at- 
mospheric contribution (jHilll . 119961 ). see Fig. [T] The neutrino generator takes 
into account the interaction of the neutrinos with the Antarctic ice. Then, 
the tracks of the neutrino- induced muons are s i mulat ed in the muon propa- 
gation Monte Carlo mrnc (IChirkin and Rhodd . |2004| ). Finally, the detector 
simulation AMASIM (iHundertmarkl . Il999[ ) is used to simulate the emitted 
Cherenkov light and to emulate the hardware behaviour. 

The same event selection applied to the experimental data set are used on 
the simulated data to identify the energy-sensitive variables for the analysis. 
The simulation data and the 10% experimental data set are compared to 
verify the consistency of the event selection. This is essentially the approach 
used in a previous AMANDA analysis of the same upgo ing neutrino data 



set to obtain limits to the extraterrestrial neutrino flux ( lAchterberg et al. 
2007ah . 



After having optimised the analysis on the 10% subset of the data, the full 
data set is used to determine the energy spectrum of the detected up-going 
neutrino-induced events. These events are expected to pile up at energies 
below 2 TeV. Givei i the high estimated puri ty of the investigated data set in 
terms of neutrinos ( lAchterberg et al.l . l2007al ). its energy spectrum is expected 
to be consistent with an atmospheric spectrum. While the theoretical pre- 
dictions of the shape of this spectrum only show minor deviations between 
2 TeV and 20 TeV, at energies between 20 TeV and 200 TeV deviations 
due to the unknown contribution of prompt neutrinos - in this case a slight 
flattening of the spectrum - are possible. Neutrinos of extraterrestrial ori- 
gin are expected to significantly flatten the spectrum at even higher energies. 



A schematic view of the analysis chain is shown in Fig. [T] Six energy- 
dependent observables are used as input variables for a neural net in order 
to produce a combined, optimised energy variable (Section E]). The com- 
bined variable is taken together with one further observable plus one of the 
variables entering the neural net, the latter still containing a component or- 
thogonal to the NN output. With the three partially independent variables, 
the neutrino energy spectrum can be determined using regularised unfold- 
ing (Section [H]). The resulting spectrum contains the effective number of 
events. As consequence of the unfolding, the probability density function of 
each event is distributed to several bins of the energy spectrum leading to 
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broken event numbers in the spectrum. By relating the simulated neutrino 
energy spectrum to the effective number of events obtained by unfolding the 
simulated flux, the energy spectrum of the effective number of events of the 
measured data is finally normalised to obtain the absolute neutrino energy 
spectrum. 

5. Selection of energy dependent variables 

Due to the limited acceptance and finite resolution of the detector, the re- 
construction of the neutrino energy spectrum can be obtained only with 
unfolding methods. For this analysis, th e meth o d of regula rised unfolding is 



used according to the RUN algorithm by iBlobell (119851 . Il996l ) . This algorithm 
allows up to three energy dependent variables as input for the unfolding pro- 
cedure. Therefore, a set of three variables correlated to the neutrino energy 
is selected. 

In a first step, the following seven observable^ which show the best correla- 
tion with the generated neutrino energy in simulations, are selected from all 
variables well described by the simulation. 

• NcHi- The number of OMs (channels) having detected exactly one 
photon as signal during the event. Due to the stochastic energy losses 
of the muons, the number of emitted Cherenkov photons and hit OMs 
increases with the muon energy. 

• NcH'- The number of OMs (channels) having detected one or more 
photons as signal, which increases with increasing muon energy. 

• Nhits'- The total number of signal photo electrons within the event, 
which increases with increasing muon energy. Each OM can contribute 
to an event by the detection of one ore more photoelectrons counted in 
this variable. 

• tmean'- The average photon arrival time, i.e. the sum of all recorded 
photon arrival times relative to the trigger time, divided by the total 



^Six observables are used as first input and an additional one is used at a later stage 
of the analysis. 
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number of hit optical modules. The higher the neutrino or muon energy, 
the more the mean time is shifted to late arrivals. 

• t JIMS'- The root mean square of the arrival time distribution of the 
photons, which grows with the number of late photons generated i.e. 
in secondary energy losses. 

• Qmean- The sum of all measured charges, in units of photoelectrons, 
divided by the number of hit OMs. This is equivalent to the mean 
number of recorded Cherenkov photons and is correlated to the energy. 

• Qrms- The root mean square of the charge distribution of the photons 
in each OM, which grows with the maximal number of photons from 
secondary energy losses. 

Some of these variables are correlated and their number exceeds the maxi- 
mum of three allowed as input for RUN. As a consequence of a multitude of 
tests, in a second step, the first six variables are combined through a neural 
network (NN) to give one energy dependent variable. The NN output is then 
used together with log^Ncn) and log{QRMs) as the three inputs to the RUN 
algorithm. 

The neural network is a standard back-propagation Multi-layer Perceptron 
(MLP) with two hidden layers used in a 6-6-3-1 feed-forward architecture 
(Fig.©. 

The complete simulation chain is used to generate a training set of muon 
events, with energies uniformly distributed in equidistant logarithmic energy 
bins between 500 GeV and 5 PeV and with trajectories uniformly distributed 
throughout a cylindrical volume with 400 m radius around the detector cen- 
ter. The NN is trained and tested with muon data sets, each containing 
100,000 events. The simulated events are then reconstructed and processed 
in the same way as the experimental data. The energy resolution of the NN 
output is estimated with four test sets of mono-energetic muons, generated 
in the same way as the training set and with muon energies of 1 TeV, 10 
TeV, 100 TeV, and 1 PeV. The resulting performance of the neural net and 
the corresponding resolutions are shown in Fig. |3l The neural net output can 
be fitted with a Gaussian distribution around the logarithm of the expected 
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Input layer 




Figure 2: Topology of the NN used for the data analysis. 



energy value. The results for the parameters of the fitted Gaussians are given 
in Table [H 



log,o(E^/GeV) 


logio(^mcan/GeV) 


logio(^/GeV) 


3.0 


3.03 


0.42 


4.0 


3.92 


0.58 


5.0 


4.99 


0.51 


6.0 


5.86 


0.48 



Table 1: Resolution of the NN output. Here, i?^, is the true input energy, i^mcan is the 
mean energy of the output and a is the standard deviation from -Bmcan- All values are 
given in logarithmic units. 



Figure |4] shows how the NN output correlates to the neutrino energy, and 
justifies its use as input for the spectrum unfolding. The comparison of 
the NN output for simulated and experimental data (Fig. [5]) shows good 
agreement. The agreement in slope depends on the neutrino energy spectrum 
chosen for the simulation. The apparently somewh at steeper dec r ease o f the 
NN output predicted for a simulation according to iHonda et al.l (120041 ) will 
be consistently visible also in the comparison of the final unfolding result 
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1 TeV 100 TeV 

lOTeV 1 PeV 




log(E^ /GeV) 



Figure 3: Output of NN for the analysis of monocnergetic muons. The muons are simulated 
with fixed energies of 1, 10, 100 and 1000 TcV. 

(Fig. [9]) with different ffux predictions. 

6. Determination of the energy spectrum 

Following the notation introduced in section |3l we form the energy sensi- 
tive variable y as the vector combination of the neutral network output, the 
logarithm of the number of channels fired and the logarithm of Qrms- The 
probability distribution of this vector variable is ^(y), leading to the need to 
solve for $(£■) in the equation p{y) = J A{y, E) ■ ^{E)dE + b{y). 

By binning the generated energy distribution and the recorded parameter, 
this integral equation can be transformed to a linear matrix equation 

y = A-E + b. (3) 

The vectors y, E and b represent the histograms containing the distribution 
of the observable, the sought-after energy spectrum and the distribution of 



15 



^ 7r 
> 

P 6.5- 



e 6: 
o 

~ 5.5 



5^ 
4-5 r 

3-5 r 
3^ 
2.5- 

4 



8 



log(E,/GeV) 



Figure 4: Correlation between NN output and true neutrino energy. Here, the mean value 
of the NN output with its errors as a function of the true neutrino energy is shown. From 
the NN output, the neutrino energy can be determined with a standard deviation of about 
0.5 order of magnitude. 



the background. The kernel A contains the design matrix, describing the 
statistical detector properties. The off-diagonal terms in the kernel arise 
from the finite resolution of the energy estimators. Solving this equation 
by inversion leads to an ill posed problem because the transfer matrix A 
necessarily contains off diagonal elements much smaller than unity, which in 
turn prevents the calculation of an a priori stable solution. To stabilise the 
solution, proper assumptions about the curvature of the solution have to be 
introduced to cut off insignificant elements of the matrix A. 



In the RUN algorithm (IBlobell . Il985l . Il996l ) . the probability distributions used 
to unfold spectra on the basis of given observed parameters are parametrised 
in the form of a superposition of cubic B-splines of fourth order. The possible 
curvature of the solution is controlled by the number of degrees of freedom 
and the number of knots of the spline-superposition. If the number of degrees 
of freedom is too small it would damp significant amplitudes and smooth the 
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Figure 5: Comparison of the neutrino energy NN output for data and a Monte Carl o 
simulation using the parametrisation of atmospheric neutrinos from (jHonda et al.l . 120041 ) . 



solution too strongly; too many degrees of freedom could enforce unphysical 
wiggles in the solution. The number of knots of the spline is of only little 
influence on the result if chosen much higher than the number of degrees 
of freedom. This procedure is called regularisation and implemented in the 
unfolding algorithm RUN. 

For the given experimental situation, the problem simplifies to the deter- 
mination of an approximately linearly decreasing function if, instead of de- 
termining $(£') from log($(i?)) is calculated from \og{y). To obtain an 
optimal parameter combination, extensive simulation tests of the following 
form were carried out: first, the transfer matrix A is calculated with an 
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arbitrar}|§ neutrino energy spectrum, a specific setting of tlie smootliing reg- 
ularisation parameter and a specific binning in y and E. The used number 



of de grees of freedom is 5 and tlie number of knots 26 (see iBlobell (11985 



19961 ) for a description of the structure). With this setting of the regularisa- 
tion parameter, a possible flattening of the spectrum is visible in all tested 
cases. The unfolding result is restricted to be positive and the RUN internal 
histogram used to calculate the acceptance correction is smoothed. 

Using these settings, different simulated spectra were unfolded. In total, the 
unfolding quality is checked with 278,000 Monte Carlo data sets. Each of 
them contains the unbiased statistical equivalent of one year of AMANDA 
data, a combination of atmospheric neutrinos and added signal contributions 
with an E'"^ spectrum. For signal contributions proportional to i?"^ ranging 
up to an contribution of 10~^ GeV cm~^ s~^ sr~^, 1,000 independent unbi- 
ased Monte Carlo data sets, each containing one year of AMANDA data are 
produced. These Monte Carlo sets were checked to see if a flattening of the 
neutrino energy spectrum towards high energies would be observed if such a 
signal were present. Also, the statistical errors obtained with the algorithm 
for the given binning follow a Gaussian (or for small effective number of 
events, a Poissonian) distribution. The chosen unfolding parameters would 
smooth out spikes on top of the spectrum smaller than the bin width in the 
middle and lower side of the investigated energy range. 

An essential part of the unfolding procedure is the proper estimation and 
accounting of the uncertainties, statistical and systematic, that propagate 
through the unfolding to determine the error bars on the spectrum. In RUN 
the statistical error is calculated under the assumption that Poissonian and 
Gaussian statistics can be applied. The same analysis shows that the distri- 
bution of the unfolding results for every bin follows a Gaussian. It further 
shows that towards high energies a flattening of the spectrum is visible with 
the chosen method if an extra signal component is present. 



In addition to the statistical errors, there are several sources of systematic 
uncertainty that affect the estimation of the unfolding matrix and thus prop- 
agate through the unfolding to the error bars on the physical atmospheric 



^It is verified that there are no relevant effects for the unfolding results, if training 
spectra and true spectra do not deviate more than ±1 in the spectral index. 
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neutrino spectrum. Estimates of the uncertain ti es due to the response of 
the detector are quantified in lAchterberg et ahl (l2007bl ). Here, we discuss 
additional sources of uncertainty from the neutrino cross sections and muon 
propagation, which both influence the rate and angular distribution of the 
detected events. 

The neutrino-nucleon DIS cross-sec tion has been measu red directly at ac- 
celerators up to only ~ 350 GeV flAmsler et al.l . |2008| ). At much higher 
energies, deep inelastic scattering probes a kinematic region (high and 
low Bjorken x) where the parton distribution functions have not been di- 
rectly meas ured. The CC v N cross sections used for t his paper are calcu- 
late d as in (IGiesel et al.l . l2003t iReya and Rodigerl . 120051 ) . with details given 



in ( iGliick et al.l . Il999l ). The calculations use the QCD inspired dynami- 
cal small X pr edictions for parto n distributions according to the radiative 
parton model ( Gliick et al.l . Il998l ) and lead to the conclusion that the sys- 
tematic uncertainty in the investigated energy range is well below 10%. 
Also, a recent exercise carried out at next-to-leading order using the ZEUS 
global PDF fits has provided the neutrino cross-section from 10'^ GeV up to 
10^^ GeV with an estimated unc ertainty ranging between ±3% and ±14% 
(ICooper-Sarkar and Sarkarl . |2008| ). 



The systematic errors due to muon energy loss, ice properties and effective 
efficiency of the photomultipliers can be estimated by comparing the recon- 
structed and expected slope of the depth intensity relation of atmospheric 



muons. 



This deviation depends on a number of factors relevant for light detection, 
including the muon energy loss, Cherenkov light propagation effects and the 
effective efficiency of the Cherenkov light detection. 

The average range of a muon is approximately R= \/h ln(l + b/a- E), where 
a and b are the effective energy loss parameters for ionisation and the sum of 
the stochastic processes (bremsstrahlung, pair production and photonuclear 
interactions) respectively. They are related to the energy loss rate through 
dE/dx = a + b ■ E. 

Although the total cross section can be calculated with high precision, the 
spectral averaged energy parameter b, in the sense used here closely connected 
to the Cherenkov light produced by the stochastic energy losses, can only be 
estimated to a precision of a few percent. 
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The probability of detecting emitted Cherenkov light depends on both on 
light propagation through the ice and the effective efficiency e of the optical 
modules. 

Comparing the slopes of depth intensity relation of atmospheric muons in 
Monte Carlo and data allows us to set an upper limit on the averaged sys- 
tematic uncertainties. In addition to the factors relevant for the muon detec- 
tion, this deviation also contains effects due to the model dependent muon 
production in the atmosphere. For this analysis, the maximal deviation in 
slope was 10%. 

This uncertainty in slope transfers directly to the neutrin o flux calculation. 



Combining all the independent detector systematics from lAchterberg et al 



( l2007bl ) of (8%) with the cross section (10%) and the muon propagation 
uncertainties (10%) gives a total uncertainty in flux of 16% which is applied 
to the statistical error bars from the unfolding. 



7. Energy spectrum of atmospheric neutrinos between 2 TeV and 
200 TeV 

In a flrst step, an intermediate energy spectrum, is determined from the data 
and compared to simulation results in Fig. |6l This result physically corre- 
sponds to the convolution of the true physical atmospheric neutrino spectrum 
and the neutrino survival and de tection efficiency, i.e. effective area of the 



detector ( lAchterberg et al.l . l2007a[ ). The final neutrino spectrum will later be 
found by correcting for the effective area and observation time. The energy 
distribution of the effective number of event^ obtained by unfolding the data 
of the years 2000-200 3 and a distribution of events simulated according to 



Volkova and Zatsepiru ( 119801 ) are compared. Shown is the number of events 
averaged in the investigated zenith angular range from 100° to 180°. The er- 
ror bars include the systematic errors of 16%. For E > 2 TeV, the unfolded 
spectrum agrees with this prediction within the errors. The corresponding 
number of events for the unfolded spectrum are given in Table HI 

In order to understand and demonstrate the effects of the energy resolution, 



^As discussed in section 6, single recorded events contribute in general to several energy 
bins of the unfolded energy spectrum. The summed weights are called effective event 
numbers. 
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Figure 6: Comparison of the energy distribution of the effective number of events obtained 
by unfolding the data of the years 2000-2003 (da ta points) and an event distribution 
simulated according to Volkova and Zatsepin ( 1980l) (histogram). 



a further check is performed. The events in the final simulated data set 
are split up into nine different energy sets based on the known true energy 
of the primary neutrino. Each of these event sets are then independently 
unfolded using the same algorithm as for the full data set. The widths of the 
resulting individual spectra give an indication of the energy resolution of the 
experiment. Each spectrum is fitted to a Gaussian 

f (log(S„/GeV)) ^ A . exp {-i (ME^/O^^^MEAjj . (4) 

The results of these fits are given in Table [31 The width of the energy bins 
denote the energy range to which the unfolding and the corresponding errors 
refer. Since the energy resolution obtained with the method described is be- 
tween 0.4 and 0.5 in \og{E,y) (see Table |3]), the contents of the bins in the 
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1 / 1 — 1 / X' Y T 7" \ 

\og{EjGeV) 


number of events 


3.3 - 3.5 


124 


3.5 — 3.7 


108 


3.7-3.9 


70.3 


3.9-4.1 


41.6 


4.1-4.3 


23.6 


4.3-4.5 


11.5 


4.5-4.7 


5.96 


4.7-5.0 


3.64 


5.0-5.5 


1.26 


5.5-6.0 


0.00 



Table 2: Effective number of events obtained by the unfolding of the neutrino data of the 
years 2000 to 2003. 



final spectrum are correlated. The statistical errors obtained by regularised 
unfolding account for this fact. With the discussed method, within the in- 
vestigated energy range a mean energy resolution of 0.45 in \og{Ey/GeV) is 
reached. 

In the final step, to obtain the actual energy spectrum of the primary atmo- 
spheric neutrinos, the "at the detector" neutrino spectrum is corrected for 
the detector efficiency and neutrino survival probability. This physical energy 
spectrum of atmospheric muon and anti-muon neutrinos is presented in Fig. 
[71 Table mhsts the values for the measured neutrino spectrum dN^/dEi, ■ 
for each energy bin. In the highest-energy bin, the error bars are compatible 
with a flux equal to zero. A flt according to Eq. ([T]), with the normalisation 
of the spectrum A and the spectral index 7 as free parameters yields 

A^ = (0.022 ± 0.026) GeV-^cm^^s^Sr"^ (5) 
7 = 2.55 ±0.13. (6) 



These values are compatible with the theoretical prediction by lVolkova and Zatsepin 
( 1l980l ). Alvoikova = 0.0285 GeV"^cm-2s-^sr-i and 7|voikova = 2.69. The er- 



ror of Ai, is compatible with zero, since it represents the flux at 1 GeV, while 
measurements are performed at above 100 GeV. 
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log(E^/GeV) 


A 


log(E^/GeV) 


a 


S.o — 3.7 


0.15 


3.56 


0.48 


3.7-3.9 


0.16 


3.76 


0.46 


3.9-4.1 


0.15 


4.01 


0.43 


4.1-4.3 


0.17 


4.22 


0.41 


4.3-4.5 


0.18 


4.41 


0.43 


4.5-4.7 


0.20 


4.61 


0.46 


4.7-5.0 


0.30 


4.92 


0.50 


5.0-5.5 


0.44 


5.39 


0.49 


5.5-6.0 


0.38 


5.73 


0.44 



Table 3: Results of the Gaussian fit to simulated data. In the first two columns the energy 
bins are given, is the normalisation, the reconstructed mean energy and a the width 
of the distribution. 
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Figure 7: The unfolded energy spectrum of muon- and anti-muon neutrinos in the atmo- 
sphere, measured with AMANDA. 
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Logytj^/ Kjey ) 


^ AT / ^ T7' ZT"^ 

ai\ / ahiy ■ tiy 

[iu >ore V / s/sr/ cm j 


6.6 — i.b 


140 27 


3.5-3.7 


83^i^ 


3.7-3.9 


47^?o 


3.9-4.1 


23^^ 


4.1-4.3 


12+3 


4.3-4.5 


5.6l^i 


4.5-4.7 


7+1.0 
^- ' -1.4 


4.7-5.0 


1 1+0.5 

-"-•-■--o.s 


5.0-5.5 


n Q4+0.20 



Table 4: The unfolded energy spectrum of muon and anti-niuon neutrinos in the atmo- 
sphere, using AMANDA data from the years 2000-2003. The errors give the 68% C.L. 
interval on the unfolded flux. 



8. Discussion of the atmospheric energy spectrum 



Figure IH] compares the unfolded energy spectrum (blue dots) to previously 
measured ener gy spectra. Measur ements by the Frejus experiment are shown 



as red s quares (Daum et al. ^ l995h. T he red lines represent SuperK measure- 
ments (iGonzalez-Garcia et al.l . |2006| ). The latter result is given in the form 



of a band of possible values indicated in this plot as an upper and a lower 
line. Allowed values (90% confidence level) lie between those two lines. The 
blue lines at higher energies represent AMANDA mea surements, based o n 



the same data sample, but optimised for low energies ( lAbbasi et al.l . l2009l ). 



Again, results are presented in form of a band, lying between the upper and 
lower line. All measurements are for the sum of neutrinos and antineutrinos. 
Although the low energy AMANDA analysis is based on the same data set 
as this analysis, it is fundamentally different from the high energy analysis 
based on an extended regularised unfolding alg o rithra discussed here. The 



low energy analysis presented in ( lAbbasi et al.l . l2009l ) used the concept of 
forward-folding. A set of curves with a limited number of parameters is 
used to give an e stimate of the i nput e nergy spect r um. I n this special case. 



the prediction by iGaisser et al.l (120011 ): iBarr et al.l ( l2004l . |2009|) is used with 
varying normalisation and spectral index to determine the spectrum. This 
method is most sensitive to the median energy of the sample, which is around 
640 GeV and is therefore not optimal for investigations at high-energies. 
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These results are the first measurement of the atmospheric neutrino spectrum 
at energies up to 200 TeV. Limits to an extraterrestrial neutrino flux with a 
generic E~'^ spectrum are shown as dashed lines in Fig. [HI The Frejus limit 
is shown at energies between 10^'^ GeV and 10^ GeV. A limit derived from 



the sam e AMANDA data set used here was presented in lAchterberg et al 



(j2007al ). confirming that no significant contribution from extraterrestrial 
sources at energies between 10^'^ GeV and 10^'^ GeV can be identified at 
the current sensitivity level. 

The unfolded neutrino energy spectrum is compared to different predictions 
of the conventional neutrino flux in Fig. [9l As the measured neutrino spec- 
trum includes zenith angles in the range 100° < 9 < 180°, the predictions 
are angle-averaged for comparison and the sum of muon and anti-muon neu- 
trinos is used. The conventional atmospheric neutrino flux depends on pa- 
rameters which le ad to uncertainties in the prediction of the expected flux 
( iBarr et al.l . |2006| ) . We compare the measured result with different predic- 
tions for this flux. The analytic approximation by IVolkova and Zatsepin 



( ll980l ) is s hown as the dot-dashed line. The solid line represents the Bartol 
pr ediction (IGaisser et al.l . 1200 ll : iBarr et al.l . |2004| . |2009| ) . The flux calculcated 
by lHonda et al.l ( 120071 ) is shown as the dashed line. The measured spectrum 
is in good agreement with all three predictions for conventional neutrinos. 

In Fig. Uni the measured energy spectrum is compared to the prediction 
of the combined spectrum of conventional and prompt neut r inos. For the 
conventional spectrum, the prediction made by iHonda et al.l (120041 ) is cho- 
sen. For the prompt contribution, several different models are shown. The 
Recombination Quark Parton M odels (RQPM) is pheno menology-based and 
non-perturbative, as described by lFiorentini et al.l ( 1200 ll ). The shown QGSM 
model is half-empiric, i.e. a combination of theoretical modeling and accel- 
erator data. This model uses the Quark Glu on Str i ng M odel based on non 



pertur bative QCD calculations, presented by I Costal (120011 ) and lBugaev et al 
(Il989l). Show n is th e maximum pr ediction. Further pr edictions are given by 



Martin et al.l (120031 ). A model by lEnberg et al.l ( 120081 ) is shown in its min- 
imum and maximum configuration. Uncertainties increase towards higher 
energi e s as e laborated in Section |2l The highest prediction {QGSM opt, 
Costa diooih ) is still compatible with the error bars of the spectrum pre- 
sented here. Next-generation experiments like IceCube will have a higher 
sensitivity to a prompt component. 
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Figure 8: The unfolded spectrum from tliis analysis (blue dots) compared to other mcasure- 
ments of the atmos pheric neutrino spectrum: Frejus measurement s (red squares) are from 



( Daum et all 119951 ) . SuperK results (red lines) are presented in ([Gonzalez-Garcia et al 



2006 ). A measurement by AMANDA at low energies (blue lines) is sho wn in (lAbbasi et al 



20091 ) ■ Limits (dashed lines) are from Frejus between lO^-^ - 10"^ GeV ()Rhode et alJ . ll996l ) 
and from a high-energy a nalysis with AMANDA using the same data set as this analysis 
(jAchterberg et all l2007at) . 



9. Conclusions and Outlook 

The unfolded muon and anti-muon neutrino energy spectrum is presented 
for the energy range 2 TeV and 200 TeV, constituting the first measure- 
ment at such high energies. The spectrum is compatible with predictions of 
the conventional and prompt atmospheric neutrino spectra. The AMANDA 
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Figure 9: Comparison of the unfolded energy spectrum to different predictions of the atmo- 
spheric neutrino flux, resulting from pion and kaon deca .ys. The predictions ar e comp atible 
with the meas ured spectrum with i n the given errors : ( Volkova and Zatsepin . 198(1) fdot- 
dashed line), (jCaisser et all I2OOII: iBarr et all |2004 l2009l ) (sohd line) and (jHonda etaT 
20071) (dashed). 



detector was switched off in May 2009 but its more than 60 times larger suc- 
cessor IceCube is currently being built at the same South Pole location. As 
of February 2010, 79 strings have been deployed and completion is planned 
within a year, completing an instrumented volume of 1 km^. 

Figure [11] shows the results of this analysis together with predictions for 
extraterrestrial neutrino fluxes. Typical neutrino fluxes from e.g. Active 
Galactic Nuclei or Gamma Ray Bursts are expected to follow a spectrum 
close to E~^, which is much harder tha n both the conventional E'^'"^) 
and t he prompt 
e.g. 



E~ ) neutrino flux (iHalzen and Hooperl . 120021 : iBeckerl . 
This implies a flattening of the spectrum towards high energies 
which is much more distinct than for prompt neutrinos. IceCube has the 



2008 
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Figure 10: The measured neutri no energy spectrum compared to predictions of a combi- 
nation of convention al neutrinos (Honda et alj . 200?!) and prompt neutrinos. The different 
prom pt models are : Fiorentini et alJ ( 2001 ): Bugaev et al.l ( 19891 ) (Bugaev RQ PM, sohd 
line): iMartin et all (l2003l) (Martin GBW dot-long-dashed hue); ICostal (l200ll) {QGSM- 
opt, dotted line) and Enberg et al. ( 20081 ) {Enberg/min, Enberg/max, dot-short-dashed 
and dashed lines) 



potential to observe this flatten ing; of the spectrum, as its main s ensitivity 
hes in the range 10^ - 10^ GeV flHoshina et al.l . l2008l : iHalzenl . l2008f ) and will 
be able to measure the high-energy neutrino spectrum with higher accuracy 
and towards higher energies than AMANDA within the flrst few years of 
operation. 
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Figure 11: Measured muon and anti-muon neutrino spectrum and predictions of ex- 
traterrest rial neutrino fluxes. Neutr in os are e xpecte d from e.g. Active Galactic Nu- 
clei (e.g. Stacker and Salaraqnl (1996h : IStecker ( 2005() . (1)), Gamma Ray Bursts (e.g. 
Waxman and Bahcalll (|l997l Il999[) ." (2)) as well as from the interactions of ultr a high - 
energy cosmic rays with the cosmic microwave background (e.g. lYiiksel and Kistler 
(3)). The expecte d sensitivity of IceCube to an £^,7^ neutrino spectrum is in the range of 
the hatched area ( Hoshina et al. . 2008 : Halzen . 20081) . 
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